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Summary 

Theoretical studies of surface roughness effects in full- 
film EHL contacts are described. The analysis, using a 
flow factor modification to the Reynolds equation, was 
applied to piezoviscous-elastic line contacts. Carefully 
converged results for ensemble-averaged film shape, 
pressure distribution, and other mechanical quantities 
were obtained. Asperities elongated in the flow direction 
by a factor exceeding 2 decreased both film shape and 
pressure extrema at constant load; isotropic or transverse 
asperities increased these extrema. Changes were small, 
of order 1 percent, and the EHL spike showed no special 
sensitivity. The largest effects were displayed by traction, 
which increased by over 5 percent for isotropic or 
transverse asperities and by slightly less for longitudinal 
roughness. 

Introduction 

The conventional elastohydrodynamic model of 
concentrated lubricated contacts is based on laminar flow 
of a Newtonian fluid between smooth, elastic bounding 
surfaces. The physical basis for this model, developed 
over the last 25 years or so, has been amply reviewed in 
the recent book of Hamrock and Dowson (ref. 1), where 
details of actual calculations and applications are given. 
With this understanding extensions of the fundamental 
model to include additional effects, such as nonlaminar 
flow, non-Newtonian rheology, or thermal behavior are 
feasible. The present work describes theoretical studies of 
some of the effects of surface roughness on the lubricant 
film, where a random texture was superimposed on the 
smooth, nominal film shape representative of the kinds 
of surface finish found in engineering practice. This work 
is particularly concerned with the random aspect of 
surface asperities. Once the local film thickness h 
becomes a stochastic variable, model solutions are 
possible only in a statistical sense, and at some point an 
ensemble averaging must be performed. Previous 
attempts to develop equations capable of handling 
roughness effects have been critically reviewed by Elrod 
(refs. 2 and 3), who notes that most early treatments were 
restricted to one-dimensional roughness textures, or 
striations, parallel to which the film profile remained 


smooth. One approach to the problem of two- 
dimensional surface roughness is to derive a modified 
Reynolds equation in which h and the pressure p are 
replaced by their ensemble-averaged values. Such an 
approach is justified to some extent since the appearance 
of h in the Reynolds equation is itself only the result of 
applying a boundary condition on the Navier-Stokes 
equation at the true surface. Its validity is, however, 
limited to surface textures that do not destroy the laminar 
flow of the film, and this sets a lower limit to the 
roughness correlation lengths measured in terms of h. 
Elrod (ref. 3) suggests 5 h as a suitable value. The validity 
of the Reynolds equation also depends on roughness 
heights being small, since otherwise pressure and velocity 
begin to vary across the film. In combination these two 
conditions amount to an average slope limitation 
distinguishing between the regimes of Reynolds and 
Stokes roughness. 

Of the many ways to derive new Reynolds equations, 
those emphasizing the conservative nature of the 
lubricant flow seem most appropriate (refs. 3 and 4). 
Such flow methods recognize that although h and p are 
local fluctuating quantities, they are also correlated since 
the flow described by their products and derivatives 
fluctuates only on a global scale. Thus, if h and p are 
replaced in the Reynolds equation by their average 
values, the effects of this correlation can be adequately 
described by a set of flow factors appearing in the 
modified Reynolds equation. These flow factors 
represent various lubricant entrainment effects of the 
surface texture and, provided that flow fluctuations are 
ignored entirely, are deterministic properties of the 
surfaces and the nominal film thickness. A calculation of 
these tensor flow factors by a perturbation expansion in 
powers of 1/A, where A is the film parameter, has been 
described elsewhere (ref. 5) in a paper referred to herein 
as I. 

The final structure of the modified Reynolds equation 
was only slightly more complex than the original smooth- 
surface form: for the case considered of a line contact in 
pure rolling with incompressible lubricant, no additional 
terms were involved. The constant fluid density, how- 
ever, was replaced by the flow factors. Consequently, 
computational methods described in the recent work of 
Hamrock and Jacobson (ref. 6), referred to herein as II, 
were readily adapted to the study of roughness effects on 



the EHL of line contacts. Although the methods are 
similar, EHL calculations are notoriously sensitive to 
changes in lubricant behavior, and the roughness, acting 
as an effective compressibility dependent on h rather than 
p, has a marked effect on the convergence of the 
solution. Following a brief summary of the formalisms 
reduced from I and II, some of the peculiarities of the 
solution procedure required by the presence of the flow 
factors are discussed. 

To the lowest non vanishing order of perturbation 
theory, roughness effects are determined by just two 
parameters. The rms surface height a referenced to h 
defines the first parameter, which in the usual notation 
becomes the film parameter A = h/a. The second 
parameter describes the anisotropy of typical surface 
asperities given by y, the ratio of the two-point 
correlation lengths parallel and perpendicular to the lay 
direction of the surface texture. For conditions typical of 
an EHL contact in the piezoviscous-elastic regime with a 
well-developed pressure spike near the outlet, the film 
shape was calculated as a function of y and A ot0 , where 
the double subscript indicates that a is normalized to the 
value of the minimum film thickness h m> 0 for the contact 
operating with a smooth surface: A m 0 =/j m 0 /<r. Other 
mechanical quantities associated with the contact have 
been computed, and typical results for the traction 
coefficient are presented. 

Symbols 

This symbol list includes definitions of some of the 
variables needed later in the formalism. The configur- 
ation to be described (fig. 1) consists of two rough, elastic 
cylinders separated by a Newtonian lubricant film and 
rotating about axes along the y direction (perpendicular 
to the page). 


Rough surface 

_ .. Nominal surface 

*" u bl 



Figure 1.— Variables used in describing fluid/solid boundary in 
lubrication of moving rough surfaces. 


b semiwidth of Hertzian contact, R V8H7rr, m 
E Young’s modulus of elasticity, Pa 

E' effective plane strain elastic modulus for two 
contacting solids, 



/ force due to shear stress on surface, N/m 
G dimensionless materials parameter, 1/Qoo 

H dimensionless nominal film thickness, h*/R 

h film thickness, m 

K dimensionless volume flow per unit length, k/uR 
k volume flow per unit length, m 2 /s 
P dimensionless nominal pressure, p*/E' 

p film pressure, Pa 

p m ax Hertzian maximum pressure, Iw/irb, Pa 

Q dimensionless isoviscous (reduced) pressure 

(eq. (14)) 

Q * dimensionless asymptotic reduced pressure 

R reduced radius of roller pair, (r - 1 + *) _ ',m 

r radius of roller, m 

S nominal dimensionless separation of undeformed 

roller pair touching at origin, x 2 /2 R 2 

t time, s 

U dimensionless entrainment velocity, riou/E'R 

u ai surface velocity component i of roller a, m/s 

«,• mean (entrainment) velocity component /, 

(u bi + u ai )/2, m/s 

Vi slip velocity component /, u bi -u ai , m/s 
W dimensionless load per unit length, w/E'R 
w load (in z direction) per unit length, N/m 

w x x component of force due to normal stress on 
surface, N/m 

X dimensionless coordinate, x/b 

x,y,z Cartesian coordinates, m 

Z Roelands viscosity-pressure exponent (eq. (16)) 

z mean height, (z&+z fl )/2, m 

z a height of surface a, m 

a viscosity-pressure index, (In rj)/p, (Pa) -1 

y pressure parameter in Roelands model (eq. (16)), 
Pa; roughness anisotropy index 

A dimensionless elastic deformation of roller pair 

subject to film pressure (eq. (12)) 

5 a random roughness height of surface a, z a -z' a , m 
6y Kronecker delta, 2x2 unit matrix 

f relative contribution of pressure gradient to flow 
(eq. (7)) 
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7} lubricant viscosity, Pa s 

rj dimensionless viscosity, 77/77 0 
rjQ viscosity at ambient pressure, Pa s 

rfi viscosity parameter in Roelands model (eq. (16)), 
Pa s 

A dimensionless film parameter, h*/a 

fx coefficient of rolling friction; traction coefficient 

v Poisson’s ratio 

a rms value of 5 b - 8 a , equal to ( 0 % + a%) 1/2 if a,b are 
uncorrelated, m 

<fPj pressure flow factor 

<p]j shear flow factor 

Subscripts: 
a roller a 

b roller b 

ij two-dimensional Cartesian vector suffixes; (x,y) 

denoted by value (1,2); implied summation on 
repeated indices 

m absolute minimum film thickness 

0 smooth surface 

Operators: 

* ensemble average (expectation) operator for 

stochastic quantity 

di,3j two-dimensional gradient operators, m - 1 

Formalism 

Flow Factor Method 

In this section the modified Reynolds equation derived 
by the flow factor method is briefly reviewed and then 
incorporated into a full computational scheme for the 
elastohydrodynamic lubrication of line contacts. For 
one-dimensional flow the formalism of I becomes 
relatively simple. 

The starting point for the flow factor method is to 
recast the flow vector kj in terms of the average pressure 
and film thickness p* and h* , which replace their true 
fluctuating values p and h. To achieve this, the effects of 
roughness must be included explicitly. For incompress- 
ible laminar flow we have adopted the form chosen by 
Patir and Cheng (ref. 4): 

*'“-w (1) 

The shear flow factor <p?j accounts for flow produced by 
roughness in the presence of slip even when the pressure 
gradient vanishes. Additional entrainment due to 


roughness under pure rolling conditions is included in the 
pressure flow factor <fP. Fluctuations in flow are not fully 
incorporated into the <p factors, but it is central to the 
method to assume that these fluctuations are negligible 
compared with those of p and h. With this assumption 
the continuity conditions d,& ( = - dh/dt can be applied to 
ki instead of k t , yielding the result 

a /(V Tif d jP*) = -Vj^i(z&u+ (2) 

where djZ is the mean gradient of the bounding surfaces 
with respect to the x,y plane. This is the modified 
Reynolds equation for surfaces in translational motion as 
given in I, where it was used to find the flow factors 
themselves in terms of the film parameter A and the 
anisotropy index 7. 

The flow factors are locally deterministic, and a 
second-order perturbation calculation using Gaussian 
forms both for the height distribution and the two-point 
autocorrelation function yielded the result 



with all other components zero. These factors depend on 
the surface height distribution only through the rms value 
a. Similarly, for given 7, the form of the correlation 
function has only a weak influence on the result. In this 
representation the x and y directions coincide with the 
roughness axes determined by the surface lay, from 
which the general case is readily obtained by coordinate 
rotation. 

The forms given by equation (3) are compared in I with 
those computed directly from an ensemble of generated 
rough surfaces (refs. 4 and 7) with generally good 
agreement (fig. 2) in the physically realistic range of A. 
This is the full-film regime, where A >3. Asperity 
contact, neglected in I, becomes important for smaller 
values, and the partial lubrication regime is generally 
taken to be 1 <A<3. By modifying the film shape near 
asperity contacts to allow for nonoverlap of the surfaces, 
the work of Bush and Hughes (ref. 8) extends the range 
of validity of the flow factor approach to A values in this 
range. 
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For rotational motion about fixed axes the form of the 
new Reynolds equation given by equation (2) requires a 
slight change since the term dh/dt in the continuity 
equation is now zero. The corresponding form of 
equation (2) then becomes 


where the subscripts have been dropped from k\ and U\. 
Applying continuity to equation (5) leads to the Reynolds 
equation appropriate to this case: 


dx \ 11 12ij dx ) 



( 6 ) 


This equation can be compared with its smooth-surface 
counterpart, equation (1) of II, the starting point of the 
full EHL calculations for line contact. For incompressible 
flow the density cancels from the equation and is replaced 
by the pressure flow factor in the Poiseuille term only. 

Equation (5) can be cast in nondimensional form by 
substituting the nondimensional variables defined in the 
symbol list, with the result 

K=(\-m (7) 

where 


r= 


J 


7T H 2 yf* dP 
2W2AU ^ dX 


( 8 ) 


The r term is the relative pressure (or Poiseuille) contri- 
bution to flow and is directly proportional to the (1,1) 
component of the flow factor. 

Similar nondijnensionalizing of equation (6) leads to 



which shows that film pressure dependency on //involves 
explicitly the two standard nondimensional groups, 
entrainment speed U and applied load per unit length W. 
Equation (9) must be solved subject to the conditions that 
P begins at ambient (effectively zero) at the inlet and that 
the Reynolds condition P = dP/dX = 0 holds at the outlet , 
the cavitation boundary. 

Nominal film thickness H is given by the sum 


d i (rfj 72^ d j P*) = “i d ' h * ~ \ ° v j d Mj W 

For pure rolling of cylinders about axes in the y 
direction, side leakage (y flow) averages to zero provided 
that the lay is either parallel or perpendicular to the rota- 
tion axis. Choosing the latter for illustrative purposes, we 
find that the flow vector of equation (1) reduces to the 
single component 


H=H 0 + S(X) + A(X) (10) 

where S(X) is the undeformed shape of the gap, A(X) is 
the elastic deformation, and H 0 is a constant equal to the 
difference between the central film thickness and the 
central elastic deformation. In the present case the 
contact has been loaded sufficiently heavily to yield 
negative values for // 0 . In the usual parabolic approxi- 
mation for circular cylinders 


, (hy J,*dp* 

'“--ST*-* +A " 


( 5 ) 


4 w 

S(X)= X 2 

7 r 


di) 
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The elastic displacement of the cylindrical boundary 
loaded by the pressure distribution P is given within the 
usual linear approximation by 

A(X) = - -J— f Pin \X-X' | dX' (12) 

7T v 7 T J 

Equations (11) and (12) together show that the depen- 
dence of H on P involves only the dimensionless load. 

Because of the importance of the pressure (and in 
general temperature)-dependent viscosity rj in the EHL 
solution, particular care must be exercised in choosing a 
model. Since r\ appears only in the denominator for dP, 
the viscosity can always be eliminated from the Reynolds 
equation by introducing the isoviscous (reduced) pressure 
Q defined by the Weibull transformation (ref. 9): 

dP 

dQ= — (13) 

V 

dP' 

<2(P)= (14) 

J i»(P') 

The viscosity model then determines the relationship of Q 
and P. For the present isothermal treatment of EHL we 
have suppressed the temperature dependence of Q and rj. 
Typically, Q reaches a finite upper limit Qca (sometimes 
known as Pi v . a s) an d rises to within a few percent of Q x 
when P reaches three or four times Q <*,. Most EHL 
contacts operate with Q close to this asymptotic value, 
which thus becomes another important EHL variable, 
conventionally taken as G = l/Q m , the materials 
parameter. 

Taking the Barus viscosity model as an example, we 
have r) = e a P, where a is the pressure-viscosity index. 
Equation (14) then gives 


fractional difference between Q and Q™ may be as little 
as 10 -20 . 

To avoid some of these difficulties, we have substituted 
the more gently varying Roelands viscosity model (ref. 
10) in place of the Barus exponential, which in any case is 
unrealistic at such high pressures. From Roelands’ 
extensive analysis we have 

i=exp([l-(l + ^) Z ]ln?lJ (16) 

where surprisingly the parameters 17 j and y are almost 
universal constants. Fluids are thus characterized by their 
individual Z and rj 0 values. For Roelands equation the 
transformation between P and Q is not particularly 
simple in either direction, and it proves more convenient 
to solve equation (9) for P directly. Thus, the reduced 
pressure no longer enters the computation explicitly. Its 
asymptotic value Q m is, however, needed to determine G, 
for which tabulations in reference 10 were used. 

Once self-consistent distributions P(X) and H(X) 
satisfying Reynolds equation (9) and the elasticity 
equation (12) have been obtained, various other 
quantities of interest can be calculated, such as the 
mechanical force components acting on the rollers and 
the friction (traction) coefficient between them. For 
example, the shear force per unit length / on cylinder a 
lies approximately in the x direction and can be written 

'•-1 (■£),.** (n > 

Integrating the Navier-Stokes equation for dp/dx across 
the film shows the integrand to be equal to - (h/2) dp/dx 
in the pure rolling condition. Since this is independent of 
z, equation (17) also gives //, and we have 


<2(P) = 


1 -e~ aEP 
uE' 


(15) 


fa~fb 5 h-f-dx 

t ij (tJC 


(18) 


from which it follows that G~aE' . 

In II the Barus fluid was used and the isoviscous 
Reynolds equation solved for Q. Equation (15) then 
yielded P, as needed for computation of the elastic 
displacement contribution to H , thus closing the iterative 
PzzH cycle. One of the biggest problems encountered in 
this approach is handling the constraint Q<Qoo when 
solving the isoviscous Reynolds equation. Pressures near 
the center of the contact may easily reach 1 GPa, 
whereupon e~ GP is of order 10“ 10 . In the spike the 
pressure can rise to more than twice this value so that the 


in which the explicit appearance of rf in equation (17) is 
now hidden. The normal pressure p acting on the two 
cylinders also has a net x component w x ~ w ax+ w bx> 
which can be expressed by 

Wx= \ h % dx (19) 

so that f a =f b - — w x /2. The shear and pressure forces 
thus balance, satisfying the equilibrium condition 
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fa + fb + w x = 0- The traction coefficient for either 
cylinder is given by 


M = 


1/1 


( 20 ) 


where the expectation value of the shear force is given in 
second-order perturbation theory by 




( 21 ) 


Application 

The influence of U, W, and G on EHL contacts was 
investigated in II. Here we focus on the effect of y and A 
contained in and as an illustrative application have 
chosen the smooth-surface conditions for case 2 of II. 
The pressure and film shape for this case are shown in 
figure 3. The pressure distribution displays a maximum 
near the Hertzian central maximum, followed by a 
definite minimum just before the spike, whose height 
( ~ 1 GPa) is more than double the central value. Values 
of the U, W, and G parameters for this case are given in 
table I, which also displays one set of raw mechanical 
data to yield these values. 

Equation (3) for the pressure flow factors shows a 
crossover from impeded to enhanced flow as y increases 


TABLE I.— MATERIAL AND LUBRICANT PROPERTIES 


Elastic modulus of steel rollers, E y Pa 2.00 x 10 11 

Poisson’s ratio for steel, v 0.3 

Inlet viscosity of paraffinic lubricant, r/ 0 , Pa s 4.1 1 x 10 ~ 2 

Constants in Roelands model: 

r?i,Pas 6.31X10- 5 

7 , Pa 1.96x108 

Roelands viscosity-pressure exponent, Z 0.67 

Dimensionless asymptotic reduced pressure, Q 2 X 10“ 4 

Effective roller radius, /?, m 1.11 xl0~ 2 

Specific loading of rollers, w, N/m 5X10 4 

Rolling velocity, u y m/s 5.94 x 10“ 1 


Hertzian semiwidth, b,m 8.02 X 10“ 5 

Hertzian maximum pressure, /? max , Pa 3.97 x 10 8 

Dimensionless EHL parameters: 

U lxlO-n 

W 2.05x10-5 

G 5 X 10 3 


through the critical value 2. Asperities, which on average 
are twice as long parallel to the flow as transverse to it, 
have no effect on mean flow. As representative of 
asperities on either side of this crossover, we have 
examined the isotropic case, 7 = 1, and the case where 
asperities are longer in the flow direction by a factor 
7 = 3. Early results showed that some significant 
departures from the smooth-surface case were present for 
unexpectedly large A ot0 values, so the range was 
extended to about 10. At the lower end, although no 
longer in the physically meaningful regime, A m 0 was 
taken as low as 1 to exaggerate small trends found at 
larger values. 
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Figure 3.— Pressure distribution and film shape for an EHL contact with operating conditions U= 10~ n , IF =2.05 x 10 -5 , and G = 5 x 10 3 . 
Elliptical Hertzian pressure distribution for same total load superimposed for comparison. (From ref. 6.) 




Computational Method 

The method is similar to that described in II and has 
been programmed according to the flowchart of figure 4. 
The equations are discretized on a uniform grid scaled 
to the Hertzian width and containing 660 points, with 
the inlet boundary condition set a distance 4 b upstream 
of the contact center. The cavitation boundary floats 
downstream near 3b/2. Reynolds equation (9) for P, 
containing the //-dependent flow factor in place of the 
P-dependent density, is transformed into an equation for 
$ =PH 3/2 and solved in loop 1 of figure 4 by iteration in 
finite difference form. The use of P rather than Q avoids 
problems of numerical precision arising in the transfor- 
mation described by equation (14). In substituting the 
Roelands for the Barus fluid this transformation is not 
easily written in closed form. Typically, the new estimate 
for 4> is underrelaxed; a weight factor of 0.1 was used for 
the present iterate, but in difficult cases it might be as low 
as 0.01. 

New P values for updating elastic displacements in 
loop 2 can likewise be weighted with previous values, but 
stability at this stage of the computation often permits a 
weight factor of 1. Elastic displacements are then 
computed by using the latest P by quadrature of equation 
(12) after analytic integration through the singularity. 
Loops 1 and 2 are continued until successive pressure 
increments fall below a chosen criterion, at which point 
the total load under the pressure curve must be compared 
with the input W value. 

In loop 3a the constant Hq of equation (10) is adjusted 
to bring the load into agreement with W, but even this 
does not produce a unique solution. Although the 
Reynolds equation is an expression of flow conservation, 
its numerical implementation by marching across the 
contact always from inlet to outlet allows systematic 
deviations from constant-flow equation (7) to appear. 
Such deviations must be suppressed by entering loop 3b, 
which again adjusts Hq. Final convergence to a unique 
solution is effected by passing continually between loops 
3a and 3b until both load and flow conservation criteria 
are simultaneously satisfied. At this point calculations of 
force components, traction, etc., are carried out and the 
program terminates. 

In the present computations convergence to the chosen 
limits proceeded more rapidly than those of II in loops 1 
and 2, and in fact it was rarely necessary to iterate more 
than once in either loop. Moreover, no tendency for 
kinks to develop near the minimum of the pressure curve 
was encountered in loops 1 and 2 at any stage. This 
convenience appears to be the combined result of using 
Roelands viscosity and pressure P as compared with 
Barus and Q in II. The worst convergence problems were 
encountered instead in loop 3a or 3b. For example, 
beginning with an initial total load W in loop 3a close to 
the input target, the pressure curve may enter an 


extensive phase in which IF increases with each iteration, 
eventually reaching some maximum value. At this point 
the film thickness is already too great, and it continues to 
rise as IF returns to target from above. Generally, W now 
becomes too low. Once this happens, loop 3b can be 
activated to improve flow conservation, which continues 
as IF passes through a minimum and approaches target 
from below. The entire cycle consisting of hundreds of 



Figure 4.— Flow chart for typical EHL line contact computation. 
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iterations may repeat several times before both loops are 
converged. 

The amplitude of the W swings can reach as much as 1 
percent of W depending on how fast H 0 is allowed to 
compensate. Since this approximates the contribution to 
W under the sharp pressure spike, it seems likely that the 
inadequacy of the numerical procedure to describe this 
spike is an important source of load oscillations. With 
present grid size the width of the spike at half-peak was 
about 1.5 spacings so that its shape, which depends on 
the viscosity model, was defined by only half a dozen or 
so points. Shifting a spike with such sharp curvature 
between adjacent nodes, with no change in its true height, 
would produce a fluctuation in estimated load of about 
0.2 percent of W. This fluctuation incidentally may be 
increased rather than reduced by using a higher order 
quadrature for W. The apparent height of the spike 
during such translation varies, of course, by a much 
larger fraction, amounting here to ~20 percent. These 
estimates agree with variations observed during a typical 
computer run, from which it is fair to assert that small 
movements of the pressure spike relative to the rollers is a 
major source of convergence difficulties. The grid at 
present is too coarse to allow an estimate of how much 
subdivision might be needed to give a reliable 
representation of the spike. 

Clearly, some aspects of the stability and convergence 
of the EHL solution still need to be better understood, 
and so long as this remains true, computational schemes 
such as those of figure 4 cannot be run fully auto- 
matically. Indeed, the computations described herein 
were implemented in the interactive mode, allowing 
control to be steered between loops by changing weight 
factors or convergence criteria. It is, however, 
particularly encouraging in this respect to find that when 
different interactive routes for the same case were 
followed, an essentially path-independent solution was 
reached. Thus, progress toward the goal of an optimized 
convergence strategy (ref. 11), if slow, has been 
nonetheless positive. 


Results and Discussion 

Flow Effects 

In guiding the computation through its final outer 
loops to convergence, the procedure was always 
terminated when the contact load agreed with the fixed 
input value to a part in 106 or better. The remaining test 
of convergence then is to examine flow variation through 
the contact, for which examples are provided in figure 5. 
For an effectively smooth surface (fig. 5(a)), where A m>0 
is larger than —10, the constancy achieved for net flow 
was comparable with that exhibited in II under actual 
smooth-surface conditions. By contrast, for the smallest 



0 80 160 240 320 400 480 560 640 

Number of nodes 


(a) Smooth boundaries. 

(b) A =0.93; 7 = 3 . 

Figure 5.— Total lubricant flow and Poiseuille relative contribution. 

A/n,0 value (0.93) in the longitudinal case (7 = 3), figure 
5(b) shows an increase in flow at the minimum film 
position due to incomplete balance between a reduced 
Couette term and an increase in the magnitude of the £ 
(or Poiseuille term). Correspondingly, for isotropic 
asperities this irregularity is inverted as a result of a 
pronounced reduction in the magnitude of £. 

According to equation (8) the Poiseuille term £ reflects 
most strongly the variation of the pressure flow factor. 
Other factors in £ are less sensitive to the roughness 
parameters. Data taken at the absolute minimum of 
figure 5 (minimum film thickness position) and 
normalized to the smooth value £ 0 are presented in figure 
6, which does indeed follow closely the pure flow factors 
of figure 2(a). 
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Figure 6.— Poiseuille flow term at minimum film thickness for rough 

surfaces (normalized to smooth-surface value -0.277). 

Film Shape Effects 

The pressure and film shape have been presented in 
figure 3 for the values of U, W, and G chosen. Since, as 
we noted earlier, roughness produces changes at the 1 
percent level, figure 3 would not be perceptibly altered 
for any of the physically reasonable values of A m>0 
considered here. On a much enlarged scale values of the 
absolute minimum film thickness H m /H m> o in the nip 
downstream from the pressure spike are shown in figure 7 
as a function of A m 0 for 7 of 1 or 3. Values are 
normalized to the smooth-surface value, which was taken 
to be the mean of the two values at the largest A m 0 
considered, differing by only 3x10-3. (This same 
normalization procedure was followed also in studying 
the variation of traction and flow.) The smooth curves 
joining the computed points have been sketched simply to 
aid visualization, and the scatter about the two lines gives 
an indication of the accuracy of the method. This is 
determined largely by the degree to which the converged 



Figure 7. — Dependence of minimum film thickness at constant load on 
film parameter for different asperity anisotropy (normalized to 
smooth-surface value // m0 = 1.99 x 10~ 5 ). 


solution is independent of the iteration path followed. As 
A m>0 becomes smaller, convergence is harder to achieve 
for any 7 value. In particular, when A <2, becomes 
negative if A 2 <3(2-7)/(l +7), which has the value 1.22 
for 7=1. The occurrence of negative flow factors near 
the minimum film thickness is an obvious indication that 
the limit of physical validity of the model has been 
surpassed. 

The film shape shown in figure 3 displays two other 
extrema: a weak maximum upstream of the spike, and a 
very shallow minimum located a few nodes downstream 
of center. Data for these secondary extrema resemble 
quantitatively those for the absolute minimum of figure 7. 
This may be compared with the results of Patir and 
Cheng (ref. 12), who performed a calculation on the inlet 
half only of a line contact to determine the effect of their 
flow factors on central film thickness. In qualitative 
agreement with figure 7, results for 7 <2 curve upward 
and vice versa, but in the isotropic case, for which 
quantitative comparison is possible, the slope of their 
curve is about 10 times greater than ours. This difference 
could well be attributed to the Grubin-like approximation 
to the inlet shape, details of which are known to have a 
large influence on the pressure buildup. 

Pressure Effects 

A discussion of pressure spike behavior has been fore- 
shadowed in the preceding section. Fractional changes in 
the central pressure maximum and in the minimum 
preceding the spike were somewhat smaller than the 
changes found in the secondary features of the film 
shape. For the most part then, P was even less sensitive to 
roughness than H. Clearly, if changes in spike height due 
to roughness also fell below 1 percent, the present 
technique could not reliably estimate the effect. All but 
two of the cases run so far were consistent with a spike of 
constant height and shape that simply moves a few tenths 
of a grid spacing upstream for longitudinal roughness or 
remains essentially stationary for the isotropic case. Only 
for the extreme A m o values did a fixed shape fail to fit. In 
the longitudinal case the height appeared to drop by 
about 5 percent for A m o~ 1, while an increase of similar 
magnitude occurred for isotropic asperities, A m Q~2. 
Although these values lie outside the physical regime for 
which the model is valid, the mathematical continuation 
is smooth, and the results are indicative of trends too 
small to observe in the physical range. Even so, the 
conclusion with respect to spike height advanced here 
depends on the assumed shape and remains subject to 
verification at finer grid sizes. 

Traction Effects 

Results of the traction coefficient calculations based on 
equation (21) have been plotted in figure 8 using the same 
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Figure 8.— -Dependence of friction coefficient on film parameter 
(normalized to smooth-surface value ^ = 5.51 x 10 -4 ). 

format as the previous figure. A significant difference 
between the two diagrams, both normalized to smooth 
values, is the general scale of the effect, which for 
traction showed changes of about 5 percent in the full- 
film range, A m>0 > 3, increasing to over 10 percent as A m0 
decreased to 2 in the partial lubrication regime. A further 
important distinction is that roughness increased traction 
for both isotropic and longitudinal asperities, and hence 
by extension for all y values, whereas the other effects 
discussed in this work each changed sign at the critical 
7 = 2 value. The large magnitude of the relative increase 
can be understood by examining equation (18), which 
shows that traction arises from the loss of x: symmetry of 
the film shape and pressure distributions as a result of 
both elastic and hydrodynamic effects — a film with a 
symmetric (e.g., Hertzian) pressure field would be 
frictionless. The influence of roughness, though small on 
H or P individually, has a large effect on the detailed 
balance of the positive and negative contributions to /*, 
given the actual unsymmetrical H and P distributions. It 
turned out in agreement with intuition that isotropic or 
transverse roughness produced the larger increase in ix. In 
the work of Bush et al. (ref. 13) on a lightly loaded point 
contact, friction increases for all y values, but in 
contrast, longitudinal textures show the larger effect. 
This difference in ordering with respect to y seems to be 
associated with the absence of side leakage in line 
contact. However, with traction in the pure rolling case 
appearing only as a result of the departure of H and P 
from perfect symmetry, the y dependence in different 
dimensions is difficult to predict. 


Conclusion and Outlook 

We have set up a formal method for incorporating the 
effects of roughness on lubricant films and have applied 
the model to the detailed computation of film shape and 
pressure distribution in an elastohydrodynamically 
lubricated (EHL) line contact. The many assumptions of 
the procedure have been discussed and include the 
following: First, roughness should be of the Reynolds 
type, which requires long wavelengths and small 
amplitudes for the asperities relative to film thickness. 
Under these conditions the flow factor modification of 
Reynolds hydrodynamics and the perturbation evalu- 
ation of these flow factors are both valid. Average flow 
has been aligned with a roughness axis, a condition 
occurring in practice in many finishing processes for 
rollers, where the lay is circumferential. By making use of 
the complete flow factor tensor, this condition can be 
relaxed to handle cases where special textures are 
imposed to achieve net transverse flow. Finally, the 
selection of operating parameters assigned the contact to 
the piezoviscous-elastic regime. Specifically, we have 
taken an incompressible, isothermal Newtonian fluid 
with Roelands viscosity entrained by two cylinders in 
pure rolling. Since the shear flow factors are also known 
from equation (3), an obvious extension of these 
computations would be to relax this last condition to 
include slip effects. 

With these basic assumptions detailed results were 
obtained for a range of film parameters A m0 
corresponding to surface roughness amplitudes from 0.1 
to 1 .0 times the smooth-surface minimum film thickness. 
Anisotropy effects were included by studying both 
isotropic asperities (7 = 1) and asperities three times 
longer on average in the flow direction than transverse to 
it (7 = 3). 

A significant conclusion of this work is the weakness of 
the influence of roughness on both film shape and 
pressure distribution. In the physically valid regime of 
A m 0 , greater than at least 2, the identifiable extrema of 
both changed fractionally by 1 percent or less. For 
longitudinal asperities both film thickness and pressure 
were reduced, with opposite changes for isotropic and 
transverse roughness. Although the two 7 values used in 
this study served to indicate the general magnitude of the 
effect, it will be worthwhile to examine the 7 dependence 
more fully, especially since there exist so many other 
treatments with which to compare the limiting one- 
dimensional roughness cases, where 7 is either 0 or 00 . 
From the pressure flow factors it seems clear that the 
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magnitude of roughness effects should increase as either 
of these extrema is approached. 

Although it might be anticipated that the height, at 
least, of the pressure spike would be sensitive to 
roughness details, the calculations again did not bear this 
out. Because of its extreme narrowness the grid size 
defined the spike with too few points, and some 
assumption about its shape was necessary. A fixed shape 
was consistent with the data, in which case the height was 
invariant with roughness to within the accuracy of the 
computation. For A m0 values at the extreme low end of 
the range a slight trend did emerge, suggesting that 
longitudinal asperities reduce the height and vice versa. 
This follows the same pattern established by the other 
extrema of P. Progress in understanding this behavior 
rests upon introducing a finer grid at least in the 
neighborhood of the spike. 

In contrast with the insensitivity of averaged film shape 
and pressure to roughness, the behavior of the traction 
coefficient manifested increases greater by almost an 
order of magnitude in the same A m> o range. This may well 
be the most important practical consequence of the 
model developed in the present work. It suggests that 
surface roughness provides further variables, in addition 
to the usual fluid rheology parameters, that could 
usefully be manipulated to control traction in an optimal 
manner. 

The present work was confined to the full-film 
lubrication regime and only hinted at effects that can be 
found in mixed lubrication. An evaluation of flow 
factors in this regime by Bush and Hughes shows that the 
largest effect of asperity contact occurs in the shear flow 
factor, suggesting that roughness effects should be 
sought in sliding situations. In addition, extension into 
the mixed regime requires a load compliance model for 
the asperities. No problems, in principle, block progress 
along these lines, but difficulties associated with 
numerical stability and convergence of such extended 
EHL computations will no doubt demand careful 
handling. The results will be of particular interest in the 
case of traction studies. 
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